Cite this Article
Guo Ping, Pan Yi-Kun, Li Long-Long, Tang Bin. Molecular dynamics simulation of decomposition and thermal conductivity of methane hydrate in porous media. Chinese Physics B, 2017, 26(7): 073101  
Permissions
Molecular dynamics simulation of decomposition and thermal conductivity of methane hydrate in porous media
Guo Ping1, †, Pan Yi-Kun1, Li Long-Long2, Tang Bin2
State Key Laboratory of Oil and Gas Reservoir Geology and Exploitation, Southwest Petroleum University, Chengdu 610500, China
School of Sciences, Southwest Petroleum University, Chengdu 610500, China

 

† Corresponding author. E-mail: guopingswpi@vip.sina.com

Abstract

The hydrate has characteristics of low thermal conductivity and temperature sensitivity. To further analysis the mechanism of thermal conductivity and provide method for the exploitation, transportation and utilization of hydrate, the effect of decomposition and thermal conductivity of methane hydrate in porous media has been studied by using the molecular dynamics simulation. In this study, the simulation is carried out under the condition of temperature 253.15 K–273.15 K and pressure 1 MPa. The results show that the thermal conductivity of methane hydrate increases with the increase of temperature and has a faster growth near freezing. With the addition of porous media, the thermal conductivity of the methane hydrate improves significantly. The methane hydrate-porous media system also has the characteristics of vitreous body. With the decrease of the pore size of the porous media, thermal conductivity of the system increases gradually at the same temperature. It can be ascertained that the porous media of different pore sizes have strengthened the role of the thermal conductivity of hydrates.

PACS: 31.15.at;44.30.+v;65.40.-b
Keyword:methane hydrate;thermal conductivity;molecular dynamics simulations;porous media
1. Introduction

Gas hydrates are inclusion compounds where hydrogen-bonded water molecules encage appropriate guest molecules (such as CH4, C2H6, C3H6, CO2, N2, and H2) under suitable pressure and temperature conditions.[1] Hydrates do not exist alone in nature, but are always accompanied by other materials such as quartz sand, ice, etc. The interest of gas hydrates research began with the discovery that the formation of hydrate could plug pipelines during the transmission of natural gas.[2] Natural gas hydrate (NGH) is a kind of abundant, widely distributed, and clean energy, which will play an important role in the future energy resources. At present, the research of natural gas hydrate including mining safety, effectiveness, and economy are still in the stage of experiment.[3–6] Gas hydrates are nonstoichiometric crystalline compounds, which belong to a group of solids called clathrate hydrates. It is known that type I, II, and H are the major crystal structures of gas hydrates, depending principally on the molecular size of guest molecules and the formation conditions.[1,7]

As we all know, the stability of hydrate changes equally sensitive to temperature, therefore the research of thermal conductivity to hydrate is very necessary. On one hand, it can provide the property physical data for gas hydrate exploitation. On the other hand, it can evaluate the quantity of hydrate resources by using the thermal conductivity of hydrate layer.[8,9] To further study the storage, transportation, and utilization of hydrate, it is completely essential to research the thermal conductivity of gas hydrate. For instance, in the design of gas hydrate storage process, the thermal conductivity is an important physical parameter. Because of the low thermal conductivity of hydrate, it is expected to enhance the thermal conductivity when using hydrate. Scientists have already carried out a lot of theoretical and experimental research on thermal performance of hydrate.

Due to the limitation of experimental methods, molecular dynamics (MD) simulation has become an important alternative method in recent years.[10–12] In the previous studies, the thermal conductivities of nanometer materials have been investigated by using MD simulation.[13,14] However, there are few MD reports related to the decomposition and thermal conductivity methane hydrate in porous media. In the paper, we performed MD simulation to study the effect of porous media on decomposition and thermal conductivity of methane hydrate. The paper is organized as follows. In Section 2, the structure models and computational methods are presented. The results concerning the thermal conductivity of pure hydrate, hydrate in porous media and hydrate in porous media with different pore sizes are given and discussed in Section 3. Finally, some conclusions are drawn from the present study in Section 4.

2. Structure models and computational methods
2.1. Simulation system and simulation method

The initial positions of the oxygen atoms in hydrate lattice were obtained from x-ray diffraction experiment,[15] and the hydrogen atoms were added in a random manner which was consistent with Bernal–Fowler rule. Figure 1 shows a model which consists of 368 water molecules and 64 methane molecules in a simulation box. As a reference system, a 2× 2× 2 unit cell replica of the cubic sI hydrate (a = 11.62 Å) was prepared with methane molecules placed at the small and large cages. The 2× 2× 3 unit cell silica ZSM-11 zeolite porous media models whose size were respectively 0.35 nm, 0.8389 nm, and 1.2053 nm were also taken into consideration, as shown in Fig. 2. The model denoted as SiO2(ZSM-11 zeolite) in the text is shown in Fig. 3. The layers are repeated due to periodic boundary conditions.

Fig. 1. Methane hydrates simulation system. Color scheme: white (H), red (O), gray (C).
Fig. 2. Different pore size of porous media models used in this paper to simulate system (a) for the pore size is 0.35 nm, (b) for the pore size is 0.8389 nm and (c) for the pore size is 1.2053 nm. Color scheme: red (O), yellow (Si).
Fig. 3. Final structure of hydrate + porous media simulation system (a) for silica (0.35 nm) + hydrate, (b) for silica (0.8389 nm) + hydrate and (c) for silica (1.2053 nm) + hydrate.

The rigid, three-center SPC/E potential model proposed by Berendsen et al.[16] was used for water-water interactions, in which each water had equilibrium O–H bond length of 1.0 Å and H–O–H bond angle of 109.5° (as listed in Table 2). The SPC/E potential is expressed as

Table 2.

Intermolecular potential parameters for the water models.[22–24]

.

Standard combination rules were used to calculate the Lennard–Jones parameters for the unlike pairs of atoms. The equations are

where rmn is the distance between molecules m and n, and is the interaction potential, qm refers to the charge on molecules m, ε0 is the permittivity of free apace, εmn and σmn are the Lennard–Jones parameters. Methane molecules are treated as spherical particles, and Lennard–Jones potentials are used for methane–methane and water–methane interactions.[17] Long-range electrostatic interactions are handled by the Ewald summation method, and nonelectrostatic interactions are truncated at a distance of 11.0 Å. The equations of motion are integrated by using the Leapfrog algorithm. The standard Lorentz–Berthelot combining rule was adopted for the methane–water nonelectrostatic interactions. The Lennard–Jones parameters of water and methane molecules used in this work are listed in tables 1 and 2. All molecular dynamics simulations are performed with periodic boundary condition by using LAMMPS molecular dynamics software package.[18]

Table 1.

Lennard–Jones parameters.[19–21]

.
2.2. Thermal conductivity calculation method of molecular dynamics simulation

Heat energy can be transmitted through solids via electrical carriers (electrons or holes), lattice waves (phonons), electromagnetic waves, spin waves, or other excitations. In metals electrical carriers carry the majority of the heat, while in insulators lattice waves are the dominant heat transporter. Normally, the total thermal conductivity λ can be written as a sum of all the components representing various excitations

where α denotes an excitation.

Thermal conductivity is defined as

where is the heat flow rate (or heat flux) vector across a unit cross section perpendicular to and is the absolute temperature. For the kinetic formulation of thermal conduction in gases, let us assume that c is the heat capacity of each particle and n is the concentration of the particles. In the presence of a temperature gradient , for a particle to travel with velocity , its energy must change at a rate of

The average distance a particle travels before being scattered is , where τ is the relaxation time. The average total heat flow rate per unit area summing over all particles is therefore

The Green–Kubo formulas relate the ensemble average of the auto-correlation of the heat flux to the thermal conductivity kappa

Combining Eqs. (5)–(8), we have

where V is the system volume, and is the heat flux.

This hydrate sample is equilibrated by using canonical ensemble (NVT) and micro-canonical ensemble (NVE) molecular dynamics methods with periodic boundary conditions at pressure 1 MPa and temperature 253.15 K–273.15 K. The result shows that the system is stable. The total simulation time is 1 ns with a time step of 1 fs, in which 500 ps is used for equilibration.

3. Results and discussion

In order to give the general features of the system configuration, snapshots of the configuration of methane hydrate dissociation process and schematic diagram of thermal conductivity in porous media at T = 273.15 K are shown in Figs. 4 and 5, respectively.

Fig. 4. Snapshots of the configuration of methane hydrate dissociation process in porous media at T = 273.15 K.
Fig. 5. Schematic diagram of thermal conductivity of hydrate in SiO2 porous media at T = 273.15 K.

In Fig. 4, with the decomposition of hydrate, methane molecules gradually accumulate on the surface of SiO2 porous media. In Fig. 5, the kinetic energy is exchanged between one or more particles in a hot layer and a cold layer at set intervals, which imposes an energy flux through the system and reduces the concentration and the temperature, and thus plays a role in the decomposition of hydrate.

3.1. The thermal conductivity of pure hydrate

The molecular structure of hydrate is very similar to that of ice, but its thermal conductivity (about 0.49 W m−1 K−1) is much smaller than that of ice (2.23 W m−1 K−1), and is closer to that of water (0.605 W m−1 K−1).[25,26] Although having a crystalline structure, the hydrate shows the thermal characteristics of vitreous body.[27] Generally the thermal conductivity of hydrate is extremely low. We have used molecular dynamics to simulate the thermal conductivity of pure hydrate, and compared the result with the actual value in literatures. The results show that the simulation value is 0.51 W m−1 K−1 and the actual value is 0.49 W m−1 K−1, which demonstrates that the model we use in this paper is reasonable. The thermal conductivity of methane hydrate at different temperatures is given in Table 3.

Table 3.

Thermal conductivity of methane hydrate at different temperatures.

.

As shown in Fig. 6, the thermal conductivity of methane hydrate increases with the increase of temperature. At the same time, the thermal conductivity increases rapidly with the increase of temperature near the freezing point, that is, the decomposition process is about to begin. It demonstrates that the lower thermal conductivity is the reason why gas hydrate exists stably. This phenomenon is also shown in Fig. 7. Figure 7 shows the change of methane hydrate thermal conductivity with time at T = 263.15 K, 268.15 K, 273.15 K. Within 1 ns, the thermal conductivity are increased by 312.04%, 125.34%, and 38.74% at 263.15 K, 268.15 K, and 273.15 K, respectively. This is consistent with the trend in literature.[28] It is well known that the formation of hydrate is accompanied by the decomposition of hydrate. Because the thermal conductivity of ice is higher than that of hydrate, the increase of the ice content in the decomposition process of pure hydrate will increase its thermal conductivity. When the temperature is 263.15 K, the thermal conductivity of hydrate increases first, then remains constant, increases again, and then remains unchanged. This phenomenon can be explained that when the decomposition rate is larger than the synthesis rate, the content of ice increases rapidly and the thermal conductivity increases. When the synthesis rate is equal to the decomposition rate, the content of ice is constant, and the thermal conductivity remains constant. The thermal conductivity remains unchanged until the hydrate is completely decomposed. Moreover, the lower the decomposition temperature is, the more easily it leads to icing. So the amount of methane hydrate decomposition at 263.15 K is much more than that of at 273.15 K in 1 ns. At a certain temperature and pressure condition, the thermal conductivity of gas hydrate increases with time.

Fig. 6. Thermal conductivity of methane hydrate at different temperatures.
Fig. 7. (color online) Variation of the methane hydrate thermal conductivity with time at T = 263.15 K, 268.15 K, 273.15 K, respectively.
3.2. The thermal conductivity of hydrate in porous media

The gas hydrate does not exist alone in nature. It is usually accompanied by other materials, such as quartz sand porous media, ice, etc. So it is important to study the thermal conductivity of hydrate in the composites.

Figure 8 compares the thermal conductivity between methane hydrate-SiO2 porous media system (in this section we consider it as a system) and pure hydrate at different temperatures. It can be seen that the thermal conductivity of methane hydrate is greatly increased by the addition of SiO2 porous media. In the case of 263.15 K, the thermal conductivity of the pure hydrate is 0.55 W m−1 K−1, while the thermal conductivity of system is 0.90 W m−1 K−1, which is increased by 72.73% compared to that of pure hydrate. This is due to the addition of SiO2 porous media will change the basic structure of hydrate. In addition, the SiO2 porous media has high thermal conductivity, it can significantly improve the heat transfer ability of the whole system and enhance energy transfer in gas hydrate, so as to avoid the icing phenomenon appeared in the late period of hydrate decomposition effectively and increase the temperature of the system, thus increasing the thermal conductivity of hydrate. For pure compound, its thermal conductivity increases with the increase of temperature, which is different from that of non-metallic crystal ordinary, shows a characteristics of vitreous body. Figure 9 shows that the thermal conductivity of methane hydrate-porous media system increases with time at a certain temperature and pressure. Hence, it can be concluded that from Figs. 8 and 9, the addition of SiO2 porous media improves the thermal conductivity of methane hydrate, but does not change the thermal conductivity of the vitreous body.

Fig. 8. (color online) The thermal conductivity of methane hydrate in porous media and the thermal conductivity of pure hydrate.
Fig. 9. (color online) Variation of methane hydrate in porous media thermal conductivity with time at T = 263.15 K, 268.15 K, 273.15 K, respectively.
3.3. Thermal conductivity of hydrate in porous media of different pore sizes

Due to the high thermal conductivity properties of SiO2 porous media, the porous media with different pore sizes were introduced into hydrate system, to strengthen the heat conduction performance of hydrate.

As shown in Fig. 10, with the decrease of the pore size of SiO2 porous media, the thermal conductivity of the system increases gradually at the same temperature. When the SiO2 porous media pore size is small, the change trend rate of thermal conductivity with temperature increases more obvious. The small size effect of porous media improves the SiO2 surface. Meanwhile, the hydrate has the characteristics of metastability. Therefore, there is a micro effect between the SiO2 porous media and the hydrate, which enhances the micro-energy transfer. From the view of hydrate, the micro effect on thermal conductivity is mainly reflected in the phonon mean free path length on.

Fig. 10. (color online) Thermal conductivities of methane hydrate in porous media of different pore sizes.
4. Conclusions

By using molecular dynamics simulation, the decomposition of methane hydrates in porous media and the influence of methane hydrates decomposition on heat conduction performance were studied. The main conclusions can be drawn as follows.

Firstly, the thermal conductivity of hydrate-SiO2 porous media system increases with the increase of temperature. At a certain temperature, the smaller the pore size, the larger enforcement impact for porous media on the thermal conductivity increase of hydrate. The main reasons for this phenomenon are the following two points. To begin with, porous media changed the basic structure of hydrate. The thermal conductivity of the porous media is several times larger than the thermal conductivity of the hydrate. With the porous media added, the energy transfer has been significantly enhanced. The small size effect of porous media improves the SiO2 surface. Meanwhile, the hydrate has the characteristics of metastability. Hence, the existing micro effect between the porous media and hydrate enhances the role of micro-energy transfer. From the view of hydrate, micro effect on thermal conductivity is mainly reflected in the phonon mean free path length on.

Secondly, the addition of porous media does not change the thermal conductivity of the vitreous body for hydrate. The thermal conductivity of methane hydrate-SiO2 porous media system also increases with the increase of temperature.

Thirdly, in the vicinity of decomposition, the thermal conductivity of hydrate increases rapidly with the increase of temperature. The low thermal conductivity of hydrate is one of the important reasons for its existence in nature. The rapidly increase of thermal conductivity means that the decomposition process is about to begin.

Acknowledgment

The authors thank the Applied Physics Research Group of Southwest Petroleum University for computational resources.

Reference
[1] Bagherzadeh S A Alavi S Ripmeester J A Englezos P 2013 Fluid Phase Equilib 358 114
[2] Max M D 2003 Natural Gas Hydrate in Oceanic and Permafrost Environments Dordrecht Kluwer Academic Publishers 9 16
[3] Holder G D Angert P 1981 Intersociety Energy Conversion Engineering Conference August 9, 1981 Atlanta, USA 810812
[4] Holder G D Kamath V A Godbole S P 1984 Annu. Rev. Energy. 9 427
[5] Kvenvolden K A 1988 Chem. Geol. 71 41
[6] Max M D Lowrie A 1996 J. Pet. Geol. 19 41
[7] Englezos P 1993 Ind. Eng. Chem. Res. 32 1251
[8] Kumar P Turner D Sloan E D 2004 J. Geophys. Res. 109 241
[9] Xu W Ruppel C 1999 J. Geophys. Res. Solid Earth 104 5081
[10] Ju Y Y Zhang Q M Gong Z Z Ji G F 2013 Chin. Phys. 22 083101
[11] Lin W Q Xu B Chen L Zhou F Chen J L 2016 Acta Phys. Sin. 65 133102 in Chinese
[12] Xu B Lin W Q Wang X G Zeng S W Zhou G Q Chen J L 2017 Chin. Phys. 26 033103
[13] Hu G J Cao B Y 2014 Chin. Phys. 23 096501
[14] Wu T Y Lai W S Fu B Q 2013 Chin. Phys. 22 076601
[15] McMullan R K Jeffrey G A 1965 J. Chem. Phys. 42 2725
[16] Berendsen H J C Grigera J R Straatsma T P 1987 J. Phys. Chem. 91 6269
[17] Jorgensen W L Madura J D Swenson C J 1984 J. Am. Chem. Soc. 106 6638
[18] Plimpton S 1995 J. Comp. Phys. 117 1
[19] Tse J S Klein M L McDonald I R 1983 J. Phys. Chem. 87 4198
[20] Abascal J L F Sanz E Fernández R G Vega C 2005 J. Chem. Phys. 122 234511
[21] Lopes P E M Murashov V Tazi M Demchuk E MacKerell A D 2006 J. Phys. Chem. 110 2782
[22] Rodger P M 1991 AIChE J. 37 1511
[23] Førrisdahl O K Kwamme B Haymet A D 1996 J. Mol. Phys. 89 819
[24] Hirai S Okazaki K Tabe Y Kawamura K 1997 Energy Convers Manage 38 S301
[25] Handa Y P Cook J G 1987 J. Phys. Chem. 91 6327
[26] Waite W F Pinkston J Kirby S H 2002 Preliminary Laboratory Thermal Conductivity Measurements in Pure Methane Hydrate and Methane Hydrate-sediment Mixtures: a Progress Report (Yokohama: Proceedings of the Fourth International Conference on Gas Hydrate) 728 733
[27] Desmedt A Bedouret L Pefoute E Pouvreau M Say-Liang-Fat S Alvarez M 2012 Eur. Phys. J. Special Topics 213 103
[28] Li D L Liang D Q Fan S S Peng H 2010 J. Nat. Gas. Chem. 19 229